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Abstract 

We consider a disordered d-wave superconductor in two dimensions. Re- 
cently, we have shown in an exact calculation that for a lattice model with 
a Lorentzian distributed random chemical potential the quasiparticle density 
of states at the Fermi level is nonzero. As the exact result holds only for the 
special choice of the Lorentzian, we employ different methods to show that 
for a large class of distributions, including the Gaussian distribution, one can 
establish a nonzero lower bound for the Fermi level density of states. The 
fact that the tails of the distributions are unimportant in deriving the lower 

bound shows that the exact result obtained before is generic. 
PACS numbers: 74.25-q, 74.25.Bt, 74.62.Dh 
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I. INTRODUCTION 



Considerable evidence for d-wave superconductivity in the high-temperature cuprate super- 
conductors has led to interest in studying the effect of disorder on d-wave paired systems. 
Unlike s-wave superconductors, where Anderson's theorem [jl| predicts negligible effect of 
nonmagnetic impurities on thermodynamic properties, simple defects are expected to be 
pairbreaking in superconductors with gap nodes, and are in fact generally thought to in- 
duce finite density of quasiparticle states A^(0) at the Fermi level. As in disordered normal 
metals, one might expect properties of such systems to depend strongly on dimensionality. 
In fact, Nersesyan et al. (NTW) have shown |^ that the usual t-matrix approximation for 
impurity scattering, which is exact in the dilute limit in 3D, breaks down for a strictly 2D 
d-wave superconductor. By mapping the problem onto a continuum model of Dirac fermions 
in a random gauge field, subsequently solved by bosonization methods, NTW claimed that 
the density of states of such systems must go to zero at the Fermi level as a power law, 
N{E) oc E°'. Later it was realized that for a realistic d-wave SC with 4 nodes on the Fermi 
surface their result might not be applicable 0. Although the real materials in question are 
quasi-2D, it is of considerable importance to establish the effect of disorder in the strictly 2D 
case because the existence of a 2D-3D crossover at low energies could invalidate the standard 
picture of low-temperature thermodynamics in a d-wave superconductor developed under 
the assumption of a finite residual density of states A^(0). 

Recently @], we have shown that for a lattice model of a disordered d-wave superconductor 
(SC) in two spatial dimensions, one can obtain an exact result for the density of states (DOS) 
N[E), provided that the disorder is modeled by a Lorentzian distribution of the chemical 
potential. The result was a finite DOS at the Fermi level, N{0)/No oc 7log4Ao/7 with No 
the normal DOS at the Fermi level, Aq the maximum value of the superconducting order 
parameter over a circular Fermi surface and 7 the width of the Lorentzian distribution. We 
also quoted rigorous lower bounds for A^(0) for a large class of disorder distributions which 
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we obtained using methods developed in a different context. One might worry that our re- 
suh for a Lorentzian distribution, while simple to obtain and exact for all energies, could be 
nongeneric, in the sense that a perturbation series based on a Lorentzian distribution cannot 
be defined due to the divergence of all moments. The proof of our lower bounds for the DOS 
in the case of more general disorder distributions then acquires a special importance. 

In this paper we therefore present in some detail the derivation of the nonzero lower bound 
for the DOS at the Fermi level which in a different context was first given in Ref. We 
stress that since our results are lower bounds, no arguments about the dependence of the 
DOS on disorder strength can be made. It suffices for our purposes to show that a lower 
bound exists, and that its existence does not depend on the specifics of the tails of the 
distribution i.e. power law decay, exponential decay or compact support of the distribution 
will all give a nonzero lower bound for the DOS. 

The paper is organized as follows: First, we formulate the problem and give a general outline 
of the proof. Second, we show preliminary calculations which will be used in the proof. We 
then derive the nonzero lower bound for a certain class of Hamiltonians. Finally, we show 
that the Hamiltonian of interest belongs to this class. We conclude with final remarks about 
cases where the method fails to give a nonzero lower bound (e.g. s-wave superconductors). 

II. FORMULATION OF THE PROBLEM AND OUTLINE OF THE PROOF 

The problem is defined by the Bogliubov- de Gennes Hamiltonian 

i/ = -(V' + /i)a3 + AV, (1) 

H being the random chemical potential, distributed according to a symmetric distribution 
function P{^r) are the Pauli matrices, (Tq is the 2x2 identity matrix). The kinetic energy 
operator -V^ is taken to act as V^\&(r) = \E'(r + 2ei) + \E'(r - 2ei) + \E'(r + 2e2) + \&(r - 2e2) 
on a function \Ef(r) of the sites r of a 2D square lattice A spanned by the unit vectors ei and 
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62- Note that this function involves displacements of two lattice sites rather than one, as 
would be the case in the simplest tight-binding representation of the lattice kinetic energy. 
For a system of fermions in the thermodynamic limit, the bare kinetic energy will then 
have a band representation quite similar to the usual tight-binding form, with no particular 
distinguishing features near the Fermi level. In the above definition of the kinetic energy, we 
have taken the hopping matrix element as our unit of energy. The bilocal lattice operator 
A'^ = Ar,r' is taken be A'^^'(r) = A[^(r + d) + - d) ± ^(r + 62) ± ^(r - 62)]. The 
Matsubara Green function G{iE) = {iEao — H)^^ determines the density of states in the 
usual way, namely by 

N{E) = —ImTr2{GrAiE ^E + ie)) (2) 

where the trace Tr2 refers to the 2x2 structure of the Hamiltonian, corresponding to quasi- 
particles and quasiholes of the superconductor. (. . .) denotes the disorder average, which 
consists of integration over the disorder variable (the chemical potential) at every site of the 
lattice, with a measure given by P{fir)dfir- 

To derive a lower bound of the average DOS we first write the Green function as 

Gm = '^'^""^^ [{H - rEa,)-^ - {H + rEa,)-') (3) 

where the new Hamiltonian is H = HDa^ = — (V^ + ^)DaQ — iADa2 has been introduced 
for formal reasons (cf. Ref. @] and below). The matrix D is diagonal with matrix elements 
Dr,r' = {—^Y^~^''^Sr,r' ('"i,''"2 are the two components of the 2D r-vector). Eq. (H) holds 
for any distribution of randomness in the Hamiltonian, before averaging. We are interested 
in the DOS at the Fermi level, i.e. at zero energy. This means that after the analytic 
continuation iE E + ie we will set E = 0. Consequently, iE — >• ie (with positive e), and 
the local DOS at the Fermi level at lattice site r reads 

Nr{0) = z}LTr2lmGrr{te) = ^Tr2[{H - ^eao);/ - {H + teao);,'] = -^Tr^iH^ + eVo);/. 

(4) 
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From the analytic properties of G it follows that is nonnegative (either positive or zero). 

The average DOS (l/|A|)(X]reA ^r) can be estimated from below using the method worked 
out in Ref . . The central idea of the proof is to divide the lattice A into finite lattice blocks 
{Sj}. We then evaluate the average DOS on these lattice blocks and also the contribution 
from the interactions of the lattice blocks. On the lattice blocks a "coarse graining" method 
will be used by relating the disorder integration over all other sites to one at the "center" 
of S. On this "center" site the range of integration of the random variable Zr = ^rD will 
be restrained to a finite interval [—a, a]. The cutoff a eliminates the contribution of the 
tails of the distribution. Since A^^ is nonnegative, the tail contribution can only add to the 
result obtained by integrating over [—a, a]. Thus, if we are able to find a nonzero average 
DOS by integrating only over the finite interval [—a, a] we have obtained a nonzero DOS 
without relying on tail contributions. This explains why distributions with power law tails 
(e.g. Lorentzian distribution) lead qualitatively to the same results as, for example, the 
Gaussian distribution or distributions with exponential decay. 

The proof rests on an identity (Eq. (14) in the next section) that is intimately connected 
with the fact that the local DOS N,. (before averaging) is nonnegative. It also relies on the 
infimum of the disorder distribution in the restricted range [—a, a] being finite. This puts 
some limits on the applicability of the proof to compact distributions, but it always holds for 
unbounded distributions like the Gaussian or the Lorentzian. The result can be summarized 
by the statement: For any finite suhregion S of the lattice A with boundary dS , defined by 
the lattice sites of A\S which are connected to S by the matrix elements of H , there exists 
a distribution dependent positive constant Ps, related to a restricted disorder distribution on 
S , with 

T^{j:Nr)>Ps{l'\dS\/\S\). (5) 

Since the block size 15*1 grows faster than the size of the boundary \dS\ the right hand side 
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is positive above a certain block size. 



It should be noted that the method of this paper will not give a nonzero lower bound for 
the DOS for every Hamiltonian. We will determine the conditions for the lower bound to 
be nonzero and show that the Hamiltonian of interest (d-wave SC) fulfills these conditions. 
We will also show that for an isotropic s-wave SC with a local order parameter the method 
will only yield a (trivial) vanishing bound for the Fermi level DOS. 



III. PRELIMINARY CONSIDERATIONS 



As a first step we evaluate the integral Nrdzr. For this purpose the identity 

-1 



Is{H + IsAIsVh = + IsAIs ^ (6) 



s 



is useful. Is is the projector onto the region S, and {.■.)s^ = 15)^^/5 is the inverse 

on the region 5*. The proof of this identity is given in Appendix A. Choosing 5* = {r}, i.e., 
just a single lattice site, we note that the 2x2 matrix ({H + ZrCJQ + iecro)^^)^^} diagonal. 
Furthermore, it is proportional to the unit matrix (Jq as a consequence of the definition of 
H, which involves only ctq and (J2. (Terms of the inverse that are proportional to (J2 are 
nonlocal and, consequently, projected out by /{r})- Therefore, we can write 

/Xr + iYr \ 

{{H + zrao + teao)-%^^ = { , (7) 

V Xr + lYrJ 

where oc e > 0. For the special choice I^r}AI{r} = —ZrO'o in Eq. (||) we obtain 

{H + ieao)„' = {{{H + z^a^ + lea^Y^}^ - z^a^X^ . (8) 
This gives for the local DOS of Eq. (^ 

Iy,,[{X,,-z,)' + Y^]-\ (9) 



TT 



The integration over Zr leads to 
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NrdZr = 1. (10) 

This result will be used below. It is nontrivial as, for example, for an s-wave SC with a local 
order parameter term Aai the equivalent of Eq. (|^ would have also off-diagonal entries. In 
fact, we have shown in |^ that the corresponding expression of the local DOS of the s-wave 
superconductor reads 

1 ie 



NriO) = — ^ ... , [(-V^ - /X - .v/A^)-i - (-V^ - /i + .VA^)-!], (11) 

i.e., it is proportional to e/VA^ + e^. This implies a vanishing DOS in the limit e — as 
long as the superconducting order parameter is nonzero. Of course, that is what is expected 
for a SC with a nonvanishing gap everywhere on the Fermi surface (Anderson's theorem) 
i- 

The local DOS Eq. (^) can be written in a differential form as 

% Q ~ ~ 

Nr = — t:; — [logdet(if — zecTo) — logdet(if + zecro)] 

4:71 OZr 

' -[logdet{l-2te{H + ieao)-% (12) 



An dZr 



which follows from the fact that the differentiation with respect to Zr picks the r, r- 
component of if ± ieao. We define a matrix A as 

A := -2ie{H + ieao)-\ (13) 

Because the DOS is nonnegative the function ilogdet(l + ^4) is a nondecreasing function of 
Zr- As a consequence of the Eq. (|T0|) we have the following integral 

t Zr = 0O 

Nrdzr = — logdet(l + A)_ =1. (14) 



oo 



47r 



-oo 



IV. LOWER BOUND OF THE AVERAGE DOS 
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A. Coarse Graining and Elimination of the Distribution Tails 



Now we consider the average local DOS summed over all lattice sites on a finite lattice block 
S, and restrict the range of integration over Zr to a finite region in order to eliminate the 
tail contributions. If we can find a nonzero lower bound for the DOS in this way, we have 
established it independently of the specific decay (e. g. power law or exponential) at large 
values of the disorder variable. First we choose a site rg G 5 for which we restrict the Zr 
integration to the interval [—a, a]. For the remaining integrations on S we define 

Zr = Zro + SZr with 6Zr G [0, 6]. (15) 

The above choice of the range of integration on S is sufficient but not necessary, i.e. different 
choices can be made as long as the range of integration is finite and certain conditions 
discussed below are satisfied. In the case considered a must be chosen large enough to 
include all singularities of the Green function. This is the case if it satisfies the inequality 

0< a- 5-4(1 + A) (16) 

as we will see below. Using the notation for this restricted averaging on 5* we have 

(Nr) > {NrYs (17) 

because A^^ is nonnegative. Then we can write with Eq.([l|) 

= l-l^fdzr,P{zr,){ n dzrP{zr))^\ogdei{l + A))^^^, (18) 

where {...)a\s refers to the (unrestricted, i.e. Zr G [—00, 00]) averaging over Zr on all lattice 
sites on A except the ones on S. 

We now bound the r.h.s. of Eq. (|18|) from below by pulling out certain infima of the 
distribution. As a first step, we pull out the infimum of the distribution on Zr^- This leads 
to 
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> inf Pizr ^ ^ 



47r 



dZr,{ n j dZrP{Zr))^\0gdet{l + A)) (19) 

h5 



Furthermore, by pulling out the infimum of the integrand of Hrgs.r^ro J'zro dzrP{zr) we 
obtain 

— iogdet(i+A)[ n 



- Pi^ro)^{ r dzr, inf -— logdet(l + A)[ [] P Pi^r)dzr 

(20) 



Pulling out the integral over z^. {tq ^ r E S) eventually yields 



> inf Pizr 



inf ( / P{Zr)dZr 
e[-a,a] ^ Jz^^ 



— ( inf r rfz,„-^logdet(l + A)\ (21) 
The disorder distribution on the block S is now taken care of by the coefficient 

, rZrQ+S , 151-1 

Ps:= inf ^P{zr^) inf (j P{zr)dzr) (22) 

•■0 



Zrf,&[-a,a\ Zrf.e\-a,a\ ^ J z. 



which multiplies the remaining disorder average over h\S. Ps is nonzero as long as the 
disorder distribution P{zr) is not vanishing in the restricted range of integration. This is 
certainly true for unbounded distributions, like a Gaussian. However, compact distributions 
with a narrow range of disorder will fail to provide a nonzero lower bound. This will be 
discussed in more detail below. 



B. General Lower Bound for the DOS 



Combining Eqs. (M), m) and (B^) we obtain 



res 



yNr) > —Ps( inf 

^ V5 - 4n \ Zr&lzo,Zo+S] 



logdet(l + A) 



log det(l + A) 



A\5 



(23) 



In the next step we isolate the lattice block S from the rest of the lattice A by sending Zr to 
±oo on the boundary dS of S. (Particles trying to occupy sites on the boundary will either 
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be trapped or repelled by an infinitely strong barrier). The boundary dS is defined by all 
sites of A which are not in S but connected with S by the matrix H, i.e, all sites r ^ S 
with \r — r'\ = 1,2 for any r' G S. (Observe that due to the definition of the Laplacian 
the "boundary" is actually two layers around the block S.) With the above definition of 
the boundary dS the matrix (1 — Iqs)H{1 — los) separates into one block matrix on S and 
another one on A\S U dS 

(1 - Ids)H{l - Ids) = IsHIs + lA\sudsHlA\suds- (24) 

Applying the identities of Appendix A, it follows that the inverse of (1 — Igs)H{l — las) 
separates into two block matrices. Consequently, lim^^^_oo(re9S') ^ = AA\ds '■= —^iei^H + 
i^cro)A\ds separates into two block matrices 

Awds = + AA\suds- (25) 

Performing the limiting process for the lower bound of the DOS we can use the fact that 
zlogdet(l + y4) is a nondecreasing function of Zr 

^ logdet(l + A) = Nr>^. (26) 



dZr 4:71 



This implies a lower bound for the r.h.s. of Eq.(PB|) if we decrease the first term in Eq.( 
by taking Zr — oo and increase the second term by taking 2:^-^00 (on the boundary of 
S). The result of this procedure is the lower bound 



A Ps{ inf 

47r \ 2^g[2o,^o+<5] 



lim logdet(l + A) 

-oo(rgf?5) 



lim logdet(l + A) _ ) (27) 



Next, we rewrite the second (negative) term by applying successively Eq. (|l^) for all r G dS. 
This yields 

lim i\ogdei{l + A) = A'K\dS\+ lim i log det(l + A). (28) 

Zr^oo{r(idS) Zr^—oo{r£dS) 



We therefore have for the expression (pTl) 
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—Ps( inf 

4:71 \ 2^g[zo,ZO+<5] 



lim logdet(l + y4) 

-oo(rG95) 



lim logdetfl + A) 

-oo{redS) 



-Ps\dS\. 



(29) 



A\(SU95) 

There is no contribution from the matrix A on A\S U dS, since this matrix part does not 
depend on = ±a. Consequently, the difference of these contributions gives zero, and we 
find a lower bound of the form 



inf 



res 



S ATT Zr€[zo,ZO+5] 



logdet(l + As 



logdet(l + A5) 



Ps\dS\. 



(30) 



The right hand side of Eq. (^) is a difference between a contribution from the block S 
(the logarithmic terms) and a boundary contribution (the \dS\ term). If the contribution 
of the block grows with its volume 151 we find for sufficiently large lattice blocks a positive 
lower bound for the r.h.s. of Eq. (|30D . We show below that this is indeed the case for the 
considered model of a d-wave SC. 



C. Lower Bound for a 2D d-wave Superconductor 

The growth of the block contribution with the volume IS*! follows from the range of the 
disorder integration on S (Eq. ([TB|) and (|T^)). To see this we define 

H' = IsHIs + ^.0^0/5. (31) 

IsHIs and H' can be diagonalized by unitary transformations. An eigenvalue \j of IsHIs 
satisfies — z^o + minAj < Aj < — z^o + maxA^. This implies for the terms in (0), where 

— — I— GE* 

=1= a + min \'j < \j < =l=a + max A^-. (32) 

An upper bound of A^-^ can be derived from the eigenvalues of Is{II + ZroDa3)Is (see Ap- 
pendix B). It yields \Xj\ < 4(1 + A) + 6, since the deterministic part of the Hamilto- 
nian — V^cts + A'^cxi has an upper bound 4(1 + A), and the random part comes from 6zr 
(0 < 6zr < 6). Thus we obtain 
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- a - 4(1 + A) < \j < -a + 4(1 + A) + (5 {zr, = a) (33) 



a 



4(1 + A) < Aj < a + 4(l + A) + 5 (^^^ = -a). (34) 



The condition for a in (W) guarantees that for ) all eigenvalues Xj are 

negative (positive). Consequently, the argument of the logarithm for any eigenvalue Xj, 
1 — 2ie/{Xj + ie), is 1 + ze (1 — ie) for the first (second) term in Eq.(P0). In order to deal 



with the branch cut of the complex logarithm we let — e ^ 27r — e for the second term in 
(0). Now we can safely let e — in both terms and obtain for 

ilogdet(l + ^5)1^ i log det{l + As)\^ (35) 

a contribution of 27r for each of the 2|S'| eigenvalues Xj, i.e. a total of 47r|S'|. From Eqs. (|l 
and (|30|) it therefore follows that the DOS is given by 



j:Nr}>Psi\S\-\dS\). (36) 

The average DOS is the sum of the local average DOS, normalized by the lattice size |A|. 
Dividing the lattice A into identical blocks S we sum over all blocks and obtain after nor- 
malization 

T^XT.Nr)>Ps{l-\dS\/\S\). (37) 

Since the lattice block size 15*1 grows faster than the size of its boundary \dS\, there is a 
finite size which gives a positive bound on the r.h.s. and therefore a positive lower bound 
on the DOS. 



Eq. ( p7D holds for our lattice model of a d-wave SC, given by the Hamiltonian Eq. (|1|) for all 
unbounded and symmetric disorder distribution that vanish at large disorder parameters Zr- 
In particular, the lower bound holds for both power law (e.g. Lorentzian) and exponential 
(e.g. Gaussian) distributions. It also holds for compact distributions of sufficient width, 
with the width being determined by the requirement that the factor Ps must be nonzero 
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when a is chosen according to the condition Eq.([l6D in order to let the DOS on S grow with 
15*1. This does not imply that narrow compact distributions will have a vanishing DOS at 
the Fermi level. However, to show the finiteness of the DOS for such distributions a more 
sophisticated method is required. 

V. CONCLUSIONS 

In conclusion, we have shown that for rather generic conditions a nonzero lower bound for 
the Fermi level density of quasiparticle states exists. The bound does not depend on the 
specifics of the "tails" of the distribution as both Lorentzian and Gaussian distributions 
yield a nonzero lower bound. This proves that our exact result for the case of Lorentzian 
disorder |^ is generic. 

This result applies to a class of Hamiltonians describing 2D superconductors with nonlocal 
order parameters, like extended s-wave, p-wave and d-wave SCs. In contrast, for a local 
isotropic s-wave SC our method will yield a vanishing lower bound, in complete agreement 
with Anderson's theorem for nonmagnetic disorder in SCs with a finite order parameter 
everywhere on the Fermi surface. It should be noted that our results imply that the self- 
consistent t-matrix approximation gives qualitatively correct physics as long as only the 
DOS at the Fermi level is concerned (i.e. for thermodynamic properties). Whether this also 
holds for the dynamic (transport) properties is an interesting question to be resolved. 
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APPENDIX A: PROJECTIONS OF THE GREEN FUNCTION 

Consider a general square matrix H defined on a lattice A. Ris a, subset of A, and Ir is the 
projector on the region R which can be written as a diagonal matrix 

r 1 iiqeR 

lR,q,q' = lR,q^qq' wlth I R,q = < ■ (Al) 

otherwise 

If the inverse of H and H + IrCIr exist then we find the identity 

{H + iRCIn)-' = - H-\l + lRCInH-%'lnCInH-\ (A2) 

where 

i...)-R' = Ir{Ir...IrVIr (A3) 
is the inverse with respect to R. From Eq. ([A2|) follows immediately 



(A4) 

and on R follows 
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Ir{H + IrCIr)-'Ir = 

by means of ([A3|) . If we choose C = Zr^ao and let — > ±00 we obtain with (A2) 

^^\imjH + IrCIr)-' = - H-\H-%'H-\ (A6) 

All matrix elements on R are zero. Therefore, we can write this expression also as a projec- 
tion onto A\R which can eventually be rewritten as the inverse on A\R 

^ (1 - Ir)H-\1 - Ir) - (1 - Ir)H-'Ir{H-YrIrH-\1 - Ir) = {H)-^{^. (A7) 

We use the above identity in the text with the choice R = dS, the boundary of the block S. 

APPENDIX B: ESTIMATION OF THE EIGENVALUES 

H and H = HDa^ are Hermitean matrices. Therefore, both matrices can be diagonalized 
by unitary transformations U and U, respectively. There are eigenvalues Xj and Xj with 

Xj = {UHU^)jj and A^- = {UHU^)jj. (Bl) 

Then we have 

A2 = ({UHU^)jj)^ = {UHWUHU^)jj = {UHDa^HD(7^U^)jj. (B2) 

Since H and Da^ commute and [Da^Y = 1? obtain for the r.h.s. 

{UH'^U^)jj<ms.-xX]. (B3) 

This estimation holds for any projection of H and if on a region S as long as the relation 
H = HDas is valid on S. We apply the above inequality in our estimation of the eigenvalues 
of the projection of H on the lattice block S. 



+ IrCIr 



R 



(A5) 
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